Code copyright the author (2022)
License:
To cite:
Code:
Inspired by https://twitter.com/DrSimEvans/status/1508409309775994892
This analysis uses the New Zealand electricity Authority Generation data.
The data contains average (mean) MW generation per half hour by source including. Where necessary we convert this to MWh per half hour (MWh = mean half-hourly MW/2)
As far as we can work out this data does not include distributed (i.e. non-grid connected) generation such as small scale wind, solar, hydro, biomass etc which is connected to the LV network. This means the EA data is likely to underestimate total generation and potentially underestimate the proportion of total generation that is renewable. It is possible that this could be fixed using embedded wind & solar generation data from the metered embedded generation data although this may still not include household level (micro)generation.
The data also includes a half-hourly carbon intensity value in g CO2/kWh sourced from https://www.carbonintensity.org.uk/.
rmdParams$olderThan <- 7
If the data we have previously downloaded is more than 7` days old, re-download.
rmdParams$dataPath <- "~/Dropbox/data/NZ_ElecAuth/processed/"
# load all the files we have using rbindlist
# probably ought to do some filtering here
filesToLoad <- list.files(rmdParams$dataPath, pattern = ".csv",
full.names = TRUE)
orig_DT <- rbindlist(lapply(filesToLoad, fread))
setnames(orig_DT, "rDateTime", "dv_dateTime")
message("Original data range from: ", min(orig_DT$dv_dateTime, na.rm = TRUE))
## Original data range from: 2020-01-01 00:15:00
message("...to: ", max(orig_DT$dv_dateTime , na.rm = TRUE))
## ...to: 2024-05-31 23:45:00
# make wide for ease of aggregation ----
nzGenMix_wide_dt <- dcast(orig_DT[!is.na(dv_dateTime)], dv_dateTime ~ Fuel_Code, value.var = "kWh", fun.aggregate = sum)
# add derived variables used later ----
addDerivedVars <- function(dt){
dt[, dv_year := lubridate::year(dv_dateTime)]
dt[, dv_date := lubridate::as_date(dv_dateTime)]
dt[, dv_month := lubridate::month(dv_dateTime)]
dt[, dv_hour := lubridate::hour(dv_dateTime)]
dt[, dv_hms := hms::as_hms(dv_dateTime)]
dt <- gridCarbon::add_season(dt,
dateVar = "dv_dateTime",
h = "S")
dt <- gridCarbon::add_peakPeriod(dt,
dateTime = "dv_dateTime")
return(dt)
}
nzGenMix_wide_dt <- addDerivedVars(nzGenMix_wide_dt)
nzGenMix_dt <- addDerivedVars(orig_DT[!is.na(dv_dateTime)])
#message("Remove incomplete years to avoid weird things in plots.")
message("Filtered data range from: ", min(nzGenMix_wide_dt$dv_dateTime))
## Filtered data range from: 2020-01-01 00:15:00
message("...to: ", max(nzGenMix_wide_dt$dv_dateTime))
## ...to: 2024-05-31 23:45:00
rmdParams$plotCap <- paste0("Data: NZ EA Generation ",
min(nzGenMix_wide_dt$dv_date), " - ",
max(nzGenMix_wide_dt$dv_date),
"\nPlot: @dataknut \nCode: https://github.com/dataknut/gridCarbon")
Note: * the data covers more years than we need * the data may contain partial years - BEWARE incomplete years in plots using annual totals or means across all months.
This looks like daily data.
# add together the % and totals we want (half-hourly)
nzGenMix_wide_dt[, dv_total := Coal + Diesel + Gas + Geo + Hydro + Solar + Wind + Wood]
nzGenMix_wide_dt[, dv_coal_gas := Coal + Gas]
nzGenMix_wide_dt[, dv_coal_gas_pc := 100 * dv_coal_gas/dv_total]
nzGenMix_wide_dt[, dv_solar_wind := Solar + Wind]
nzGenMix_wide_dt[, dv_solar_wind_pc := 100 * dv_solar_wind/dv_total]
# keep the vars we want for clarity
temp <- nzGenMix_wide_dt[, .(dv_dateTime, dv_coal_gas, dv_solar_wind,
dv_coal_gas_pc, dv_solar_wind_pc, dv_total)]
temp[, dv_date := lubridate::date(dv_dateTime)]
# aggregate to daily data for plotting
plotDT <- temp[,
.(mean_dv_solar_wind_pc = mean(dv_solar_wind_pc),
mean_dv_coal_gas_pc = mean(dv_coal_gas_pc),
total_dv_coal_gas = sum(dv_coal_gas),
total_dv_solar_wind = sum(dv_solar_wind),
total_dv_total = sum(dv_total),
nObs = .N), # to check for days with < 48 half hours
keyby = .(dv_date)
]
plotDT[, dv_year := lubridate::year(dv_date)] # for plots
plotDT[, total_dv_coal_gas_pc := 100 * total_dv_coal_gas/total_dv_total] # daily %
plotDT[, total_dv_solar_wind_pc := 100 * total_dv_solar_wind/total_dv_total]
message("Check for days with less than 48 hours - this will be truncated data due to DST breaks. We hate DST breaks")
## Check for days with less than 48 hours - this will be truncated data due to DST breaks. We hate DST breaks
table(plotDT$nObs)
##
## 46 48
## 4 1609
Figure ?? shows the mean half-hourly % generation by each type per day. This is slightly convoluted - it is the mean of the sum of the 48 daily half-hourly values. Unfold the code above for clarity.
The smoothed curves are estimated for each year. The lines terminate at the maximum value for the year. I’m still trying to decide if they tell us anything useful.
ggplot2::ggplot(plotDT[dv_year > 2011], aes(x = mean_dv_solar_wind_pc,
y = mean_dv_coal_gas_pc,
colour = as.factor(dv_year),
alpha = dv_year)) +
geom_point() +
geom_smooth() +
scale_colour_viridis_d(name = "Year") +
guides(alpha = "none") +
labs(x = "Solar & wind (mean % of half-hourly generation per day)",
y = "Coal & gas (mean % of half-hourly generation per day)",
caption = rmdParams$plotCap)
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'
Figure 4.1: Mean half-hourly % generation by each type per day
# save it
ggplot2::ggsave(filename = "meanNZrenewablesVsfossilHalfHourPC.png",
path = rmdParams$plotPath,
height = 5)
## Saving 7 x 5 in image
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'
Figure ?? shows the percentage of daily generation by type. This is less convoluted as it is the sum of generation per day for the two categories (solar + wind vs gas + coal) as a % of total daily generation.
Again the smoothed curve is estimated for each year.
ggplot2::ggplot(plotDT[dv_year > 2011], aes(x = 100 * total_dv_solar_wind_pc,
y = 100 * total_dv_coal_gas_pc,
colour = as.factor(dv_year),
alpha = dv_year)) +
geom_point() +
geom_smooth() +
scale_colour_viridis_d(name = "Year") +
guides(alpha = "none") +
labs(x = "Solar & wind (% of total daily generation)",
y = "Coal & gas (% of total daily generation)",
caption = rmdParams$plotCap)
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'
Figure 4.2: Percentage of daily generation by type
ggplot2::ggsave(filename = "dailyNZpcGenMix.png",
path = rmdParams$plotPath,
height = 5)
## Saving 7 x 5 in image
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'
Just cos we can… helpfully split into ‘peak’ and ‘off peak’ periods.
Peak period definitions:
Again the smoothed curve is estimated for each year (and demand period).
ggplot2::ggplot(nzGenMix_wide_dt[dv_year > 2011], aes(x = dv_solar_wind_pc,
y = dv_coal_gas_pc,
alpha = dv_year,
colour = as.factor(dv_year))) +
geom_point() +
facet_wrap(. ~ dv_peakPeriod) +
geom_smooth() +
scale_colour_viridis_d(name = "Year") +
guides(alpha = "none") +
labs(x = "Solar & wind (% of half-hourly generation)",
y = "Coal & gas (% of half-hourly generation)",
caption = rmdParams$plotCap)
## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'
Figure 4.3: Percentage of half-hourly generation by type
ggplot2::ggsave(filename = "halfHourlyNZ_PCgenByPeakPeriod.png",
path = rmdParams$plotPath,
height = 5)
## Saving 7 x 5 in image
## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'
Values are in kWh per half hour https://www.emi.ea.govt.nz/Wholesale/Datasets/Generation/Generation_MD
plotDT <- nzGenMix_dt[,
.(dailyMWh = sum(as.numeric(kWh/1000))), # MWh -> MWh
keyby = .(dv_year, dv_date, dv_season, Fuel_Code)]
# annual sum for cross-check
t <- plotDT[, .(sum_TWh = sum(dailyMWh/1000000)), keyby = .(Year = dv_year)]
t[, Year := as.factor(Year)]
gridCarbon::make_flexTable(t, caption = "Total TWh per year for cross-check")
Year | sum_TWh |
|---|---|
2020 | 40.6 |
2021 | 41.0 |
2022 | 40.7 |
2023 | 40.6 |
2024 | 16.9 |
# double check
Figure 5.1 shows individual trends.
p <- ggplot2::ggplot(plotDT, aes(x = dv_date, y = dailyMWh/1000, colour = Fuel_Code)) +
geom_line() +
facet_grid(Fuel_Code ~ .) +
scale_color_viridis_d(name = "Source") +
labs(x = "Date",
y = "GWh",
caption = rmdParams$plotCap)
p
Figure 5.1: Trends in daily GWh generation
plotly::ggplotly(p)
Figure 5.1: Trends in daily GWh generation
ggplot2::ggsave(filename = "dailyNZGenBySource.png",
path = rmdParams$plotPath,
height = 5)
## Saving 7 x 5 in image
Figure 5.2 stacks them
ggplot2::ggplot(plotDT, aes(x = dv_date, y = dailyMWh/1000, fill = Fuel_Code)) +
geom_col(position = "stack") +
scale_fill_viridis_d(name = "Source") +
labs(x = "Date",
y = "GWh",
caption = rmdParams$plotCap)
Figure 5.2: Trends in daily GWh generation (stacked)
ggplot2::ggsave(filename = "genNZTrendStackBySource.png",
path = rmdParams$plotPath
)
## Saving 7 x 5 in image
Defined how?
Renewable = wind + solar + hydro + geothermal + wood - but note that geothermal is not necessarily low carbon
nzGenMix_wide_dt[, RENEWABLE := Geo + Hydro + Solar + Wind + Wood]
nzGenMix_wide_dt[, RENEWABLE_perc := 100*RENEWABLE/dv_total]
make_NgesoPlots <- function(dt,
var = "RENEWABLE_perc", # defaults
lowColour = "green",
highColour = "red",
scaleLab = "Change me!",
minMax = "max"){
res <- list() # for the results
res$tile <- ggplot2::ggplot(dt, aes(x = dv_date,
y = dv_hms, fill = get(var))) +
geom_tile() +
theme(legend.position = "bottom") +
scale_x_date(date_labels="%b %Y",date_breaks ="12 month") +
theme(axis.text.x = element_text(angle = 60, hjust = 1)) +
scale_fill_continuous(name = scaleLab,
high = highColour,
low = lowColour) +
labs(x = "Date",
y = "Half-hour",
caption = rmdParams$plotCap)
# if(minMax = "max"){
# res$tile <- res$tile +
# scale_fill_continuous(name = paste0("Maximum ", scaleLab)
# )
# }
# if(minMax = "min"){
# res$tile <- res$tile +
# scale_fill_continuous(name = paste0("Half-hourly ", scaleLab)
# )
# }
# line plot
plotDT <- dt[,.(max = max(get(var)),
mean = mean(get(var)),
min = min(get(var))
),
keyby = .(dv_date)]
p <- ggplot2::ggplot(plotDT, aes(x = dv_date)) +
theme(legend.position = "bottom") +
scale_x_date(date_labels="%b %Y",date_breaks ="12 month") +
theme(axis.text.x = element_text(angle = 60, hjust = 1))
if(minMax == "max"){
yh <- max(plotDT$max)
p <- p + geom_line(aes(y = max, colour = max)) +
scale_color_continuous(name = paste0("Maximum half-hourly ", scaleLab),
high = highColour,
low = lowColour) +
geom_smooth(aes(y = max), colour = "grey") +
geom_hline(yintercept = yh, colour = highColour) +
annotate("text", x = mean(plotDT$dv_date),
y = yh*1.05, label = paste0("Maximum: ", round(yh)))
p <- p + labs(x = "Date",
y = paste0("Maximum half-hourly ", scaleLab),
caption = rmdParams$plotCap
)
}
if(minMax == "min"){
yh <- min(plotDT$min)
p <- p + geom_line(aes(y = min, colour = min)) +
scale_color_continuous(name = paste0("Minimum half-hourly ", scaleLab),
high = highColour,
low = lowColour) +
geom_smooth(aes(y = min), colour = "grey") +
geom_hline(yintercept = yh, colour = lowColour) +
annotate("text", x = mean(plotDT$dv_date),
y = yh*1.05, label = paste0("Minimum: ", round(yh)))
p <- p + labs(x = "Date",
y = paste0("Minimum half-hourly ", scaleLab),
caption = rmdParams$plotCap
)
}
res$line <- p
return(res)
}
rpc <- make_NgesoPlots(nzGenMix_wide_dt,
var = "RENEWABLE_perc",
lowColour = "red",
highColour = "green",
scaleLab = "Renewable generation (%)",
minMax = "max" # max line or min line?
)
rpc$tile
rpc$line
## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'
rpc <- make_NgesoPlots(nzGenMix_wide_dt[dv_year > 2022],
var = "RENEWABLE_perc",
lowColour = "red",
highColour = "green",
scaleLab = "Renewable generation (%)",
minMax = "max"
)
rpc$tile +
scale_x_date(date_labels="%b-%Y",date_breaks ="3 month")
## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
rpc$line +
scale_x_date(date_labels="%b-%Y",date_breaks ="3 month")
## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'
wpc <- nzGenMix_wide_dt[, .(mean_pc = mean(RENEWABLE_perc)), keyby = .(dv_hms,
Year = dv_year)]
ggplot2::ggplot(wpc, aes(x = dv_hms, alpha = Year, group = Year, y = mean_pc)) +
geom_line() +
labs(x = "Time of day",
y = "% renewable generation")
#plotDT <- gbGenMix_dt[, .()]
nzGenMix_wide_dt[, WIND_perc := 100*Wind/dv_total]
ci <- make_NgesoPlots(nzGenMix_wide_dt,
var = "WIND_perc",
lowColour = "red",
highColour = "green",
scaleLab = "Wind %",
minMax = "max"
)
ci$tile
ci$line
## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'
wpc <- nzGenMix_wide_dt[, .(mean_pc = mean(WIND_perc)), keyby = .(dv_hms,
Year = dv_year)]
ggplot2::ggplot(wpc, aes(x = dv_hms, alpha = Year, group = Year, y = mean_pc)) +
geom_line() +
labs(x = "Time of day",
y = "% wind generation")
ci <- make_NgesoPlots(nzGenMix_wide_dt[dv_year > 2022],
var = "WIND_perc",
lowColour = "red",
highColour = "green",
scaleLab = "Wind %",
minMax = "max"
)
ci$tile +
scale_x_date(date_labels="%b-%Y",date_breaks ="3 month")
## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
ci$line +
scale_x_date(date_labels="%b-%Y",date_breaks ="3 month")
## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'
#plotDT <- gbGenMix_dt[, .()]
nzGenMix_wide_dt[, SOLAR_perc := 100*Solar/dv_total]
ci <- make_NgesoPlots(nzGenMix_wide_dt,
var = "SOLAR_perc",
lowColour = "red",
highColour = "green",
scaleLab = "Solar %",
minMax = "max"
)
ci$tile
ci$line
## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'
wpc <- nzGenMix_wide_dt[, .(mean_pc = mean(SOLAR_perc)), keyby = .(dv_hms,
Year = dv_year)]
ggplot2::ggplot(wpc, aes(x = dv_hms, alpha = Year, group = Year, y = mean_pc)) +
geom_line() +
labs(x = "Time of day",
y = "% solar generation")
ci <- make_NgesoPlots(nzGenMix_wide_dt[dv_year > 2022],
var = "SOLAR_perc",
lowColour = "red",
highColour = "green",
scaleLab = "Wind %",
minMax = "max"
)
ci$tile +
scale_x_date(date_labels="%b-%Y",date_breaks ="3 month")
## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
ci$line +
scale_x_date(date_labels="%b-%Y",date_breaks ="3 month")
## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'
#plotDT <- gbGenMix_dt[, .()]
nzGenMix_wide_dt[, HYDRO_perc := 100*Hydro/dv_total]
ci <- make_NgesoPlots(nzGenMix_wide_dt,
var = "HYDRO_perc",
lowColour = "red",
highColour = "green",
scaleLab = "Solar %",
minMax = "max"
)
ci$tile
ci$line
## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'
wpc <- nzGenMix_wide_dt[, .(mean_pc = mean(HYDRO_perc)), keyby = .(dv_hms,
Year = dv_year)]
ggplot2::ggplot(wpc, aes(x = dv_hms, alpha = Year, group = Year, y = mean_pc)) +
geom_line() +
labs(x = "Time of day",
y = "% hydro generation")
ci <- make_NgesoPlots(nzGenMix_wide_dt[dv_year > 2022],
var = "HYDRO_perc",
lowColour = "red",
highColour = "green",
scaleLab = "Hydro %",
minMax = "max"
)
ci$tile +
scale_x_date(date_labels="%b-%Y",date_breaks ="3 month")
## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
ci$line +
scale_x_date(date_labels="%b-%Y",date_breaks ="3 month")
## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'
#plotDT <- gbGenMix_dt[, .()]
nzGenMix_wide_dt[, GEO_perc := 100*Geo/dv_total]
ci <- make_NgesoPlots(nzGenMix_wide_dt,
var = "GEO_perc",
lowColour = "red",
highColour = "green",
scaleLab = "Solar %",
minMax = "max"
)
ci$tile
ci$line
## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'
wpc <- nzGenMix_wide_dt[, .(mean_pc = mean(GEO_perc)), keyby = .(dv_hms,
Year = dv_year)]
ggplot2::ggplot(wpc, aes(x = dv_hms, alpha = Year, group = Year, y = mean_pc)) +
geom_line() +
labs(x = "Time of day",
y = "% geothermal generation")
ci <- make_NgesoPlots(nzGenMix_wide_dt[dv_year > 2022],
var = "GEO_perc",
lowColour = "red",
highColour = "green",
scaleLab = "Geothermal %",
minMax = "max"
)
ci$tile +
scale_x_date(date_labels="%b-%Y",date_breaks ="3 month")
## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
ci$line +
scale_x_date(date_labels="%b-%Y",date_breaks ="3 month")
## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'
XX need to create a half-hourly
CARBON_INTENSITYvariable XX Requires: * per-fuel CI
“The grid-average emission factor best reflects the carbon dioxide equivalent emissions associated with the generation of a unit of electricity purchased from the national grid in New Zealand in 2020”
“We calculate purchased electricity emission factors on a calendar-year basis and based on the average grid mix of generation types for calendar years. The emission factor accounts for the emissions from fuel combustion at thermal power stations and fugitive emissions from the generation of geothermal electricity. Thermal electricity is generated by burning fossil fuels.
The emission factor for purchased grid-average electricity does not include transmission and distribution losses. ”
“This emission factor also doesn’t reflect the real-world factors that influence the carbon intensity of the grid such as time of year, time of day and geographical area. Therefore, a grid -average emission factor may over or underestimate your organisation’s GHG emissions.”
2020 value = 0.120 kg CO2/kWh
Instead use a model derived from the GB MWh generation mix for 2023 (which may be way off but…)
Coefficients: * Intercept: 41.231171 * GAS_perc: 3.571570 * COAL_perc: 10.050410 * WIND_perc: -0.451285 * HYDRO_perc: 1.147793 # why is hydro +ve? * SOLAR_perc: -0.335122
This is (obviously) completely broken as it does not include a coefficient for Geothermal or Wood and the GB mix includes Nuclear, Interconnection and Storage. #YMMV
Need to fix :-)
# should be annual and fuel specific
# 2020 value
nzGenMix_wide_dt[, COAL_perc := Coal/dv_total]
nzGenMix_wide_dt[, GAS_perc := Gas/dv_total]
nzGenMix_wide_dt[, COAL_perc := Coal/dv_total]
nzGenMix_wide_dt[, CARBON_INTENSITY := 41.231171 +
GAS_perc * 3.571570 +
COAL_perc * 10.050410 +
WIND_perc * -0.451285 +
HYDRO_perc * 1.147793 + # why is hydro +ve?
SOLAR_perc * -0.335122]
summary(nzGenMix_wide_dt$CARBON_INTENSITY)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 66.26 100.38 106.94 106.58 113.15 136.09
Seasons:
plotDT <- nzGenMix_wide_dt[, .(mean_CI = mean(CARBON_INTENSITY)), keyby = .(dv_year, dv_date, dv_season)]
ggplot2::ggplot(plotDT[dv_year > 2011], aes(x = dv_date, y = mean_CI, colour = dv_season)) +
geom_point() +
scale_color_viridis_d(name = "Season") +
theme(legend.position = "bottom") +
labs(x = "Date",
y = "Mean g CO2e/kWh",
caption = rmdParams$plotCap)
Figure 6.1: Mean half-hourly carbon intensity per day
ggplot2::ggsave(filename = "meanCiTrendBySeason.png",
path = rmdParams$plotPath
)
## Saving 7 x 5 in image
Re-draw as a boxplot of mean daily CI by month - plotted at month start. 6.2 shows outliers nicely.
plotDT <- nzGenMix_wide_dt[, .(mean_CI = mean(CARBON_INTENSITY)),
keyby = .(dv_year,
dv_date,
dv_season)]
plotDT[, dv_month := lubridate::floor_date(dv_date, unit = "months"),]
ggplot2::ggplot(plotDT[dv_year > 2011], aes(x = dv_month,
y = mean_CI,
group = dv_month,
colour = dv_season)) +
geom_boxplot() +
scale_color_viridis_d(name = "Season") +
theme(legend.position = "bottom") +
labs(x = "Month",
y = "Mean g CO2e/kWh",
caption = rmdParams$plotCap)
Figure 6.2: Mean daily CI boxplots per month
ggplot2::ggsave(filename = "meanCiTrendByMonthSeason.png",
path = rmdParams$plotPath
)
## Saving 7 x 5 in image
Figure 6.3 shows half-hourly carbon intensity over time.
#plotDT <- nzGenMix_wide_dt[, .()]
ci <- make_NgesoPlots(nzGenMix_wide_dt,
var = "CARBON_INTENSITY",
lowColour = "green",
highColour = "brown",
scaleLab = "Carbon intensity (g CO2/MW)",
minMax = "min"
)
ci$tile
Figure 6.3: Half-hourly carbon intensity over time
ci$line
## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'
ci <- make_NgesoPlots(nzGenMix_wide_dt[dv_year > 2021],
var = "CARBON_INTENSITY",
lowColour = "green",
highColour = "brown",
scaleLab = "Carbon intensity (g CO2/MW)",
minMax = "min"
)
ci$tile +
scale_x_date(date_labels="%b-%Y",date_breaks ="3 month")
## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
ci$line +
scale_x_date(date_labels="%b-%Y",date_breaks ="3 month")
## Scale for x is already present.
## Adding another scale for x, which will replace the existing scale.
## `geom_smooth()` using method = 'loess' and formula = 'y ~ x'
Mean carbon intensity per year and season within peak period.
Figure 6.4 suggests evening peak periods still have slightly higher carbon intensity and the shape of the reduction curves differ by season although rather less by period. Interestingly the sustained reduction in carbon intensity in Summer has leveled off.
plotDT <- nzGenMix_wide_dt[, .(mean_CI = mean(CARBON_INTENSITY)), keyby = .(dv_year, dv_peakPeriod, dv_season)]
ggplot2::ggplot(plotDT, aes(x = dv_year, y = mean_CI, colour = dv_peakPeriod)) +
geom_line() +
scale_color_viridis_d(name = "Peak period") +
facet_wrap(. ~ dv_season) +
labs(x = "Year",
y = "Mean g CO2e/kWh",
caption = rmdParams$plotCap)
Figure 6.4: Mean half-hourly carbon intensity by peak period and season 2012-2022
ggplot2::ggsave(filename = "annualMeanCIByPeak.png",
path = rmdParams$plotPath
)
## Saving 7 x 5 in image
Do we see a relationship between renewables generating and peak demand? This will be mediated by the way the electricity market works.
We may find wind curtailment (not visible here) at low demand periods where nuclear can’t be shut off.
First, what is the general average shape of carbon intensity and renewable generation?
Figure 7.1 shows that although the mean half-hourly carbon intensity had fallen over time (with the effect of solar in summer particularly noticeable), the morning and evening peaks are still relatively more carbon intense as demand overtakes the available renewable supply.
plotDT <- nzGenMix_wide_dt[dv_year > 2011, .(mean_ci = mean(CARBON_INTENSITY),
mean_renewables_MW = mean(RENEWABLE),
mean_renewables_pc = mean(RENEWABLE_perc)),
keyby = .(dv_year, dv_hms, dv_season)]
ggplot2::ggplot(plotDT, aes(x = dv_hms, y = mean_ci,
alpha = dv_year,
colour = dv_year,
group = dv_year)) +
geom_line() +
scale_alpha_continuous(name="Year") +
scale_color_continuous(name = "Year",
low = "grey",
high = "#3CBAC6") + # UoS from Marine palette
theme(axis.text.x = element_text(angle = 45, hjust=1)) +
facet_wrap(. ~ dv_season) +
labs(x = "Time of day",
y = "Mean half-hourly carbon intensity",
caption = rmdParams$plotCap)
Figure 7.1: Mean half-hourly carbon intensity by year and season
ggplot2::ggplot(plotDT, aes(x = dv_hms, y = mean_renewables_MW/1000,
colour = dv_year,
alpha = dv_year,
group = dv_year)) +
geom_line() +
scale_alpha_continuous(name = "Year") +
scale_color_continuous(name = "Year",
low = "grey",
high = "#3CBAC6") + # UoS from Marine palette
facet_grid(dv_season ~ .) +
labs(x = "Time of day",
y = "Mean renewables (GW)",
caption = rmdParams$plotCap)
Figure 7.2: Mean half-hourly renewable generation by year and season
ggplot2::ggplot(plotDT, aes(x = dv_hms, y = mean_renewables_pc,
colour = dv_year,
alpha = dv_year,
group = dv_year)) +
geom_line() +
scale_alpha_continuous(name = "Year") +
scale_color_continuous(name = "Year",
low = "grey",
high = "#3CBAC6") + # UoS from Marine palette
facet_grid(dv_season ~ .) +
labs(x = "Time of day",
y = "Mean % renewables (%)",
caption = rmdParams$plotCap)
Figure 7.3: Mean half-hourly renewable generation by year and season
ggplot2::ggplot(nzGenMix_wide_dt, aes(x = RENEWABLE, y = dv_total)) +
geom_point() +
facet_wrap(. ~ dv_year) +
theme(axis.text.x = element_text(angle = 45, hjust=1))
7.4 shows the rise of both wind and solar (mid-day peak)…
Beware incomplete yearss
nzGenMix_wide_dt[, GENERATION := dv_total]
plotDT <- nzGenMix_wide_dt[, .(mean_renewables = mean(RENEWABLE),
mean_solar = mean(Solar),
mean_generation = mean(GENERATION)),
keyby = .(dv_year, dv_hms, dv_peakPeriod, dv_season)]
ggplot2::ggplot(plotDT[dv_year > 2016], aes(x = dv_hms,
colour = dv_peakPeriod,
alpha = dv_year,
group = dv_year)) +
geom_line(aes(y = mean_renewables/1000)) +
scale_alpha_continuous(name = "Year") +
scale_color_discrete(name = "Peak period") +
theme(axis.text.x = element_text(angle = 45, hjust=1)) +
facet_wrap(. ~ dv_season) +
labs(x = "Time of day",
y = "Mean renewables (GW)",
caption = rmdParams$plotCap)
Figure 7.4: Trends in mean half-hourly renewable generation by time of day and season
What would have happened if we increased solar generation in 2022 by a factor of 10? We would have needed some storage in Spring & Summer…
ggplot2::ggplot(plotDT[dv_year == 2022], aes(x = dv_hms,
colour = dv_peakPeriod,
group = dv_peakPeriod)) +
geom_point(aes(y = mean_generation/1000)) +
geom_line(aes(y = mean_solar/1000 * 10), linetype = "dotdash") +
scale_color_discrete(name = "Peak period") +
theme(axis.text.x = element_text(angle = 45, hjust=1)) +
facet_wrap(. ~ dv_season) +
labs(x = "Time of day",
y = "Mean generation (GW)",
caption = rmdParams$plotCap)
Figure 7.5: Comparing mean total generation & 10 * solar generation by time of day for 2020
YMMV on 7.6
plotDT <- nzGenMix_wide_dt[, .(mean_renewables = mean(RENEWABLE),
mean_generation = mean(GENERATION)),
keyby = .(dv_year, dv_hms, dv_peakPeriod, dv_season)]
ggplot2::ggplot(plotDT[dv_year > 2016], aes(x = mean_generation/1000 , y = mean_renewables/1000,
colour = dv_peakPeriod)) +
geom_point() +
scale_color_discrete(name = "Period") +
facet_grid(dv_season ~ dv_year) +
labs(x = "Mean total generation (GW)",
y = "Mean renewables (GW)",
caption = rmdParams$plotCap)
Figure 7.6: Mean renewables vs Mean total generation
That’s it.
You might want to look at recent academic research on this topic:
skimr::skim(nzGenMix_wide_dt)
| Name | nzGenMix_wide_dt |
| Number of rows | 77416 |
| Number of columns | 32 |
| Key | dv_dateTime |
| _______________________ | |
| Column type frequency: | |
| Date | 1 |
| difftime | 2 |
| factor | 2 |
| numeric | 26 |
| POSIXct | 1 |
| ________________________ | |
| Group variables | None |
Variable type: Date
| skim_variable | n_missing | complete_rate | min | max | median | n_unique |
|---|---|---|---|---|---|---|
| dv_date | 0 | 1 | 2020-01-01 | 2024-05-31 | 2022-03-17 | 1613 |
Variable type: difftime
| skim_variable | n_missing | complete_rate | min | max | median | n_unique |
|---|---|---|---|---|---|---|
| dv_hms | 0 | 1 | 900 secs | 85500 secs | 42300 secs | 48 |
| hms | 0 | 1 | 900 secs | 85500 secs | 42300 secs | 48 |
Variable type: factor
| skim_variable | n_missing | complete_rate | ordered | n_unique | top_counts |
|---|---|---|---|---|---|
| dv_season | 0 | 1 | FALSE | 4 | Aut: 22080, Sum: 20208, Win: 17664, Spr: 17464 |
| dv_peakPeriod | 0 | 1 | FALSE | 5 | Ear: 22582, Day: 22582, Eve: 12904, Lat: 12896 |
Variable type: numeric
| skim_variable | n_missing | complete_rate | mean | sd | p0 | p25 | p50 | p75 | p100 | hist |
|---|---|---|---|---|---|---|---|---|---|---|
| Coal | 0 | 1 | 97716.97 | 93256.93 | 0.00 | 0.00 | 79688.00 | 161903.00 | 365201.00 | ▇▅▂▂▁ |
| Diesel | 0 | 1 | 361.80 | 3439.98 | 0.00 | 0.00 | 0.00 | 0.00 | 79213.00 | ▇▁▁▁▁ |
| Gas | 0 | 1 | 254654.88 | 125511.52 | 18133.00 | 167618.50 | 245349.00 | 334624.50 | 646409.00 | ▃▇▆▂▁ |
| Geo | 0 | 1 | 424445.46 | 31472.62 | 224851.00 | 409428.00 | 430211.00 | 445384.50 | 548599.00 | ▁▁▃▇▁ |
| Hydro | 0 | 1 | 1379334.99 | 343833.61 | 492522.00 | 1110369.25 | 1396932.00 | 1626214.50 | 2363701.00 | ▂▆▇▅▁ |
| Solar | 0 | 1 | 79.75 | 749.31 | 0.00 | 0.00 | 0.00 | 0.00 | 10920.00 | ▇▁▁▁▁ |
| Wind | 0 | 1 | 152139.19 | 96910.12 | 0.00 | 70422.75 | 138378.00 | 221762.25 | 446238.00 | ▇▇▆▃▁ |
| Wood | 0 | 1 | 13058.15 | 5454.13 | 0.00 | 12047.00 | 15243.00 | 16652.25 | 19561.00 | ▂▁▁▆▇ |
| dv_year | 0 | 1 | 2021.73 | 1.29 | 2020.00 | 2021.00 | 2022.00 | 2023.00 | 2024.00 | ▇▇▇▇▃ |
| dv_month | 0 | 1 | 6.19 | 3.47 | 1.00 | 3.00 | 6.00 | 9.00 | 12.00 | ▇▆▅▅▆ |
| dv_hour | 0 | 1 | 11.50 | 6.92 | 0.00 | 5.00 | 11.00 | 17.00 | 23.00 | ▇▇▆▇▇ |
| dv_total | 0 | 1 | 2321791.20 | 401583.73 | 1324732.00 | 1983700.50 | 2353673.00 | 2600489.75 | 3626574.00 | ▂▆▇▂▁ |
| dv_coal_gas | 0 | 1 | 352371.85 | 181359.86 | 19460.00 | 208045.00 | 332673.00 | 482707.00 | 962020.00 | ▆▇▆▂▁ |
| dv_coal_gas_pc | 0 | 1 | 14.89 | 6.71 | 1.02 | 9.61 | 14.81 | 19.97 | 35.68 | ▃▇▇▃▁ |
| dv_solar_wind | 0 | 1 | 152218.95 | 96918.56 | 0.00 | 70534.75 | 138392.00 | 221797.25 | 449415.00 | ▇▇▆▃▁ |
| dv_solar_wind_pc | 0 | 1 | 6.74 | 4.47 | 0.00 | 3.08 | 5.97 | 9.74 | 24.57 | ▇▆▃▁▁ |
| RENEWABLE | 0 | 1 | 1969057.55 | 336540.82 | 1143537.00 | 1701580.50 | 1988560.00 | 2214708.75 | 3062302.00 | ▂▆▇▃▁ |
| RENEWABLE_perc | 0 | 1 | 85.10 | 6.72 | 64.32 | 80.01 | 85.19 | 90.39 | 98.98 | ▁▃▇▇▃ |
| WIND_perc | 0 | 1 | 6.74 | 4.47 | 0.00 | 3.08 | 5.97 | 9.74 | 24.57 | ▇▆▃▁▁ |
| SOLAR_perc | 0 | 1 | 0.00 | 0.03 | 0.00 | 0.00 | 0.00 | 0.00 | 0.47 | ▇▁▁▁▁ |
| HYDRO_perc | 0 | 1 | 58.90 | 7.97 | 28.86 | 53.82 | 59.03 | 64.19 | 83.17 | ▁▂▇▆▁ |
| GEO_perc | 0 | 1 | 18.88 | 3.81 | 9.56 | 16.00 | 18.05 | 21.60 | 34.20 | ▂▇▅▂▁ |
| COAL_perc | 0 | 1 | 0.04 | 0.04 | 0.00 | 0.00 | 0.03 | 0.07 | 0.21 | ▇▃▂▁▁ |
| GAS_perc | 0 | 1 | 0.11 | 0.05 | 0.01 | 0.08 | 0.11 | 0.14 | 0.25 | ▃▇▇▃▁ |
| CARBON_INTENSITY | 0 | 1 | 106.58 | 9.92 | 66.26 | 100.38 | 106.94 | 113.15 | 136.09 | ▁▂▇▇▁ |
| GENERATION | 0 | 1 | 2321791.20 | 401583.73 | 1324732.00 | 1983700.50 | 2353673.00 | 2600489.75 | 3626574.00 | ▂▆▇▂▁ |
Variable type: POSIXct
| skim_variable | n_missing | complete_rate | min | max | median | n_unique |
|---|---|---|---|---|---|---|
| dv_dateTime | 0 | 1 | 2020-01-01 00:15:00 | 2024-05-31 23:45:00 | 2022-03-17 12:00:00 | 77416 |
Packages etc: